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Abstract: We calculate the next-to-leading order (NLO) quantum chromodynamics (QCD) corrections to double 
charmonium production processes e + e~ — — > g c +h c {lP) within the non-relativistic QCD (NRQCD) fac¬ 
torization framework. We find that the corrections to r/ c +h c (lP) production are positive, while those to rjc+ipl, 2 {lD) 
are negative. Unlike the J/ip + g c case, all the corrections here are not large. Uncertainties in the renormalization 
scale, quark mass and running energy of center-of-mass are discussed, and the scale dependence of these processes is 
found to be greatly reduced with the NLO QCD corrections. 
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1 Introduction 


Once, the double charmonium inclusive(exclusive) 
processes e + e~ —> J/ip + cc(r] c ) encouraged many inves¬ 
tigations because of large discrepancies between leading 
order (LO) calculations [l|-|4[ and experimental results at 
B-factories HI . It was found that this problem may 
be somehow remedied by next-to-leading order (NLO) 
Quantum Chromodynamics (QCD) corrections J7|-| 11 ] in 
the framework of non-relativistic QCD (NRQCD)|12j . 
which indicates that the NLO corrections might be sig¬ 
nificant to LO results. In order to further understand the 
NLO properties and investigate the convergence of per¬ 
turbative expansion for heavy quarkonium production 
and decays in the NRQCD formalism, we need to study 
more heavy quarkonium production and decay processes. 

There have been great advances in recent years in the 
calculation of radiative corrections to charmonium inclu¬ 
sive and exclusive production and decays. The LO esti¬ 
mations for h c inclusive production at B-factories were 
given in Refs. ESQ. Recently, a complete NLO cal¬ 
culation for the P-wave charmonium Xcj^Pj * 1 ^, 3 M*') in¬ 
clusive production processes e + e~ -^Xcj + cc/gg/qq{q = 
u,d,s) was carried out 0- Though so far there are in¬ 
sufficient data on the double excited charmonium exclu¬ 
sive processes e + e~ -A 7 * —> H 1 + hR theoretical studies 
have already started 10 , and even the NLO correc¬ 


tions for J/ip+Xcj production at B-factories [17j,ll8j have 
been performed. 


In this work, we calculate the NLO QCD corrections 
for e + e~ -A rj c + h c (lP)/ip lt2 (lD) processes within the 
NRQCD formalism. The /i c (l 1 Pi, J PC = l +_ ) state was 
first observed by experiment E760 at Fermilab (0, and 
its C-parity was established by radiative decay h c —> ? 7 C 7 
0. The ^(3770) has been identified as ip 1 (l 3 D 1 , J PC = 
1 __ ) since its parameters are consistent with the expec- 
taions (2lj| . Note, the 0(3836) was once considered to be 
the quark model state i/> 2 (l 3 .D 2 , J PG = 2 __ ) [231, how¬ 
ever, the ■0(3823) recently observed by Belle [23J and 
BESIII 241 turns to be the ^ 2 (1 3 .D 2 ) state. For simplic¬ 
ity, we denote 0i ;2 ( 1 3 ^i, 2 ) states as ?/q, 2 throughout this 
paper. After the running of super B-factories in future, 
the processes considered in this study might be observed, 
and a precise evaluation is hence necessary. 

The rest of this paper is organized as follows: Sec¬ 
tion 2 presents our formalism and calculation method, 
numerical evaluation and some discussion of the results 
are given in Section 3, and Section 4 gives a summary 
and conclusions. 


2 Calculation Scheme Description and 
Formalism 


In our calculation the Mathematica package Fey- 
nArts 0 was applied to generate all the Feynman dia¬ 
grams and amplitudes of partons. The standard projec- 
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tion operators for charnronia may be expressed as 
1 


flop — 


4\/2 E{E + m c ) 

1 


WNl ' 


(p-m c ) </ 0A (f + 2E)(p+m c ) 

( 1 ) 


Here, n 0 corresponds to the spin-singlet charnronia with 
$) = 7 5 , while n x corresponds to spin-triplet states with 
ft = ft, the spin polarization vector, p = £ + q and 
p = — q are the momenta of the quark and antiquark 

within the charmonium, respectively. P denotes the mo¬ 
mentum of quarkonium and q is the relative momentum 
between quarks inside the quarkonium. 1 stands for the 
unit matrix in color space. 



Fig. 1. The typical Feynman diagrams of processes 
7*-^c + MhP)M, 2 (TD) in NLO QCD. 


Using the method described in Ref. ij, and apply¬ 
ing FeynCalc |27j to assist the calculation of amplitudes, 
one can readily obtain the tree-level results. For one- 
loop QCD corrections, the representative Feynman di¬ 
agrams are shown in Fig. 1. We first used FeynCalc 
to trace the Dirac matrix chains as well as the color 
matrices, and to perforin the derivative on momentum 
q. Next, the package $Apart Q was employed to re¬ 
duce the propagators of each one-loop diagram into lin¬ 
early independent ones. Then, we applied the package 
FIRE [29] to reduce all the one-loop integrals into several 
master-integrals(A 0 , D 0 .C 0 ,D 0 ). Finally, package Loop- 
Tools [3Qi was employed to evaluate the scalar master- 
integrals numerically. 

Throughout our calculation, we adopted Feynman 
gauge, and took the conventional dimensional regulariza¬ 
tion with D = 4 — 2e to regularize the divergences. The 
ultraviolet divergences are canceled by the counter terms 
and the infrared divergences in the short distance coef¬ 
ficients cancel each other, which confirms the NRQCD 
factorization for e+e - —> rj c + h c (lP)/ip li2 (^D) processes 
at NLO level. When handling the counter terms, we 
found that terms related to Z 3 , the renormalization con¬ 
stant corresponding to the gluon field, vanish in the 
end. The renormalization constant Z g , corresponding 


to the strong coupling constant a s , was computed in 
the modified-minimal-subtraction (MS) scheme, while 
Z 2 and Z, n , corresponding respectively to the quark field 
and quark mass, were in the on-shell (OS) scheme. Even¬ 
tually, the expressions for the relevant renormalization 
constants read: 


6Z™ S = 


6Z? S = 


ft) ^ s 
ft a 9 

— Of — 
4-7T 

a, 


-7.e + hi(47r) 

e uv 

1 2 47TU 2 

-1-37b + 3 In —-—f- 4 

euv e m to 2 


SZ° s = —3 a 


4-7T 


1 47r/i 2 

-7e +m- 5 - 

euv m 


( 2 ) 


3 Numerical results and analysis 

Before carrying out numerical calculation, the input 
parameters need to be fixed. The NLO running coupling 
constant 

_ 1 ft In L 

47t ft L ft 3 L 2 [ ’ 

was employed with L = ln(p 2 / Aq CD ), ft = 
(ll/3)C A -(A/3)T F n f and ft = (3A/3)C\~ AC F T F n f - 
(20/3 )C A T F rif. In numerical evaluation, the effective 
quark flavor number n.f = 4 was adopted and A q C d = 
0.297GeV fin]]. At leading order in relativistic expan¬ 
sion M r)c / hc /^ 12 = 2m c , and hence the charm quark mass 
m c = 1.7±0.2GeV was taken. The values of NRQCD ma¬ 
trix elements were evaluated from the Bunchmiiller-Tye 
potential model (32], i.e., 

N 

mr, c ™^\ R sm 2 = 0.387 GeV 3 , 

3N , 

(Oi) hc ~^\R P (0)\ 2 = 0.107 GeV 5 , 


(OiK,^^f\R'nm 2 = 0.054 GeV 7 . (4) 


Table 1. Cross sections cr(/r)(fb) of e+e - —> r/ c + 
h c (lP) at leading order and next-to-leading or¬ 
der at (i = 2m c and y's/2. The errors are induced 
by m c = 1.7±0.2 GeV. 


e+e -^r] c + h c (lP) 


cr (/r)(fb) 

LO NLO 

K-factor 

p = 2 m c 

n 970 —0.143 /D qq- 1 —0.154 
U.Z, / O_|_o.3i5 U.OO±_|_o.309 

1 1Q + 012 
i.iy_o.n 

p = yfsj2 

o.207;“+° 7 o.29i;°is? 

1 4D +011 


Table 2. Cross sections cr(/r)(fb) of e+e - —> r/ c + 
ft (ID) at leading order and next-to-leading or¬ 
der at p = 2 m c and \fs/2. The errors are induced 
by m c = 1.7±0.2 GeV. 
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e + e — >t] c 

+ ^i(l D) 


<Km)C fb) 

LO 

NLO 

K-factor 

fj, = 2 m c 

D 91 7 -0124 
U.Zl 1 _|_o .273 

0.129+8;?If 

0.596;?;??? 

H = y/s/2 

0 1fi9-° 093 

U. 1 UZi_(_ q 203 

0.144;?;?®® 

0.890;?;??? 

Table 3. Cross sections a(/i)(fb) of e + e“ —> r/ c + 

■ 02 (IP) at leading order and next-to-leading or¬ 
der at p = 2 m c and \fs/2. The errors are induced 
by m c = 1.7± 0.2 GeV. 


e + e — >r) c 

+ ^ 2(1 D) 


cr(n)(fb) 

LO 

NLO 

K-factor 

H = 2 m c 

0.869;?;“? 

0.605;?;??? 

0.6961?;??? 

[1 = yfs/2 

0.648+g;||| 

0.632;?;??? 

0.9761?;??? 


In our calculation, two typical renormalization scales 
were considered, therefore the corresponding values of 
running coupling constant are a s (/x = 2 m c ) = 0.235 
and a a (/j, = y/s/2) = 0.203. The cross sections of 
e + e“ -» rj c +h c (lP) /ip lt2 {lD) are presented in Tabs. OH 
in which the errors are induced by varying m c , and the 
K-factor is defined as ——. These results indicate that: 

G LO 

(1) For r] c + h c (lP) production, the NLO correction en¬ 

hances tree-level result, and the K-factor is bigger 
at large scale \i = y/s/ 2 ; 

( 2 ) For ’q c + i/j 1 ^(lD) productions, the NLO corrections 

are negative, and the K-factors are also bigger at 
large scale /j, = yfs/2\ 

(3) At LO, cross sections of these three processes de¬ 

crease with the scale /j. increasing from 2 m c to 
y/s/2. The NLO result of r] c + h c (lP) also de¬ 
crease with the scale, while the NLO results of 
^c+Vh, 2 ( 1 . 0 ) have inverse relation with the energy 
scale in the discussed region; 

( 4 ) The hierarchy a^ 1 < <Jh c < <Jy, 2 remains for both LO 

and NLO results at scales /j = 2m c and [7 = y/s/ 2. 

To illustrate the results in Tables m more clearly, 
we show the LO and NLO cross sections of ?; c + 
/i c (lP)/'i/'i i 2 (lTl) production with different m c versus 
renormalization scale /r in Figs. EH respectively. One 
may notice that the scale dependence of the NLO results 
is obviously depressed, and the convergence of pertuba- 
tive expansion for these processes works well, whereas 
the cross sections are quite sensitive to the quark mass. 



g/GeV 

Fig. 2. LO and NLO cross sections with different 


me versus running scale [i for e + e — > r/ c +hc(lP). 



Fig. 3. LO and NLO cross sections with different 
m c versus running scale fi for e + e _ —> ri c +tpi(lD). 



m c versus running scale p for e + e —> r/ c +ip 2 (lD). 

In Figs. Ed we present the LO and NLO cross sec¬ 
tions for e + e _ —> r) c + h c (lP)/ip 1} 2 {lD) versus center-of- 
mass energy ( E cm =y/s ), respectively, where scale p, and 
mass m c are fixed. The figures show that NLO correc¬ 
tion for rj c +h c (lP) production is positive in the displayed 
region, while that for r] c + is negative. Interest¬ 

ingly, for the r] c + 1 ^ 2(1 D) channel, the NLO correction 
has an inversion at E cm = 10.36 GeV. Going across this 
energy point, the NLO correction changes from positive 
to negative. Generally, all curves go down with the in¬ 
crease of center-of-mass energy. It is worth noting that 
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NRQCD factorization formalism for double charmonium 
only holds when y/s m c [33], and it may break down 
at the double charmonium threshold because of the color 
transfer effect [34j, whereas in drawing Figs. 03 the 
center-of-mass energy has been extended to the thresh¬ 
old region, lower than the typical B-factory energy, just 
for a schematic display. 



Fig. 5. LO and NLO cross sections versus run¬ 
ning energy of center-of-mass(-Ecm) for e + e“ —> 
Vc + h c (lP) with m c — 1.7 GeV at /x = 5.3 GeV. 



Fig. 6 . LO and NLO cross sections versus run¬ 
ning energy of center-of-mass(i? cm ) for e + e“ —» 
Vc + ipi(lD) with rric = 1.7 GeV at n = 5.3 GeV. 



Fig. 7. LO and NLO cross sections versus run¬ 
ning energy of center-of-mass(i5 cm ) for e + e“ —» 
* 7 c + ^ 2(1 D) with to 0 = 1.7 GeV at n = 5.3 GeV. 


Since 2 > M D p, i.e. it is above the open charm 
threshold, and r] c +/ lr2 (lD) cross sections are at the order 
of 1 fb, in experiment it turns out to be harder to detect 
ip lt 2 than h c through double charmonium processes at 
super B-factories. For r] c + h c (lP) production, defining 
the relative production ratio R as 

g[e+e~ —> r]c + h c \ 

a[e+e~^r] c + J/rpy {> 

and considering the Belle data of tr[e + e _ —> rj c +J/ip\ = 33 
fb [Ej], we have R = 0.010(0.009) at ft = 2m c (y / s/2). If the 
BaBar data cr[e+e _ —> rj c +J/il>\ = 17.6 fb [g] are adopted, 
we then have R = 0.019(0.017) at /z = 2m c (y/s/2). If 
we adopt the NLO results of Ref. @ for cr[e + e _ —>• 
r] c +J/ip] = 18.9 fb, then 1? = 0.018 under the same renor¬ 
malization scale /r = 2m c . In future super B-factories, 
this ratio R may stand as a benchmark for the estima¬ 
tion of the possibility of observing the h c (lP) through 
double charmonium process. 

4 Summary and Conclusions 

In this work, we have studied the double charmonium 
processes e + e _ —> rj c +h c (lP) at NLO accuracy 
under the NRQCD factorization mechanism. Cross sec¬ 
tions with varying charm quark mass m c = 1.7±0.2 GeV 
at typical renormalization scales (/r = 2m c ,y/s/2) were 
analyzed in detail. The magnitudes of cross section ver¬ 
sus energy scale /i and center-of-mass energy y/s at LO 
and NLO were evaluated. We have also estimated the rel¬ 
ative production rate of R = <r[e + e _ —> r] c +h c ]/a[e + e~ —>■ 
which might be helpful for the measurement of 
double charmonium exclusive production in future super 
B-factories. 

Through our study, we find that the NLO correc¬ 
tions for e + e _ —► rj c +h c (lP) are small, and the 
convergence of perturbative expansion works well up to 
NLO. When the scale /j, lies in the region [2m c ,y/s/2], 
the NLO correction for rj c + h c (lP) production is posi¬ 
tive, yet it is negative for r 7 c +' 0 i. 2 (1-D). The scale depen¬ 
dence has been clearly reduced when NLO corrections 
are taken into account. At y/s = 10.6 GeV, the relation 

< &h c < 0 V 2 exists at both LO and NLO, no matter 
whether /r = 2 m c or yL = y/s/2. It is worth noting that 
the cross sections of the three processes considered are 
sensitive to the charm quark mass, which hence is the 
main source of uncertainty in our results. 

Although h c (lP) has been measured through the 
h c —> J/i/> + 7 T° —> (e + e“) + 7 T° process, its branch¬ 
ing fraction is not determined, and there has been no 
signal observed in the other h c primary decay mode 
h c —> J/iI>+2tt —> (e + e _ )+27r 0 - The results in this work 
may be helpful to h c and NRQCD studies in B-factories 
in future. 
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